﻿using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace SystemAnalizer.Algorithms
{
    public class Faddeev
    {
        private double[, ,] Discret_matrix;
        private int system_size;
        private int discrete_num;
        private double[, , ,] A_matrix;
        private double[,] P_trace;
        private double[, , ,] B_matrix;

        public Faddeev(int size, int numberdsc, double[, ,] A_disc)
        {
            this.system_size = size;
            this.discrete_num = numberdsc;
            this.Discret_matrix = new double[this.discrete_num, this.system_size, this.system_size];
            this.A_matrix = new double[this.system_size + 1, this.discrete_num, this.system_size, this.system_size];
            this.B_matrix = new double[this.system_size, this.discrete_num, this.system_size, this.system_size];
            this.P_trace = new double[this.system_size, this.discrete_num];

            for (int l = 0; l < this.discrete_num; l++)
                for (int i = 0; i < this.system_size; i++)
                    for (int j = 0; j < this.system_size; j++)
                    {
                        this.Discret_matrix[l, i, j] = A_disc[l, i, j];
                        this.A_matrix[0, l, i, j] = A_disc[l, i, j];
                    }

        }

        public double[,] calculate_coef()
        {
            for (int l = 0; l < this.system_size; l++)
            {
                for (int k = 0; k < this.discrete_num; k++)
                {
                    calcktrace(l, k);
                    calck_Bmatrix(l, k);
                    calk_Amatrix(l, k);
                }
            }

            return this.P_trace;
        }

        private void calk_Amatrix(int iter_num, int discrete_num)
        {

            double currentsum = 0;
            for (int l = 0; l <= discrete_num; l++)
                for (int i = 0; i < this.system_size; i++)
                    for (int j = 0; j < this.system_size; j++)
                    {
                        currentsum = 0;
                        for (int n = 0; n < this.system_size; n++)
                            currentsum += this.A_matrix[iter_num, discrete_num - l, i, n] * this.B_matrix[iter_num, l, n, j];
                        this.A_matrix[iter_num + 1, discrete_num, i, j] += currentsum;
                    }
        }



        private void calck_Bmatrix(int iter_num, int discrete_num)
        {
            for (int i = 0; i < this.system_size; i++)
                for (int j = 0; j < this.system_size; j++)
                    this.B_matrix[iter_num, discrete_num, i, j] = this.A_matrix[iter_num, discrete_num, i, j];

            for (int i = 0; i < this.system_size; i++)
                this.B_matrix[iter_num, discrete_num, i, i] -= this.P_trace[iter_num, discrete_num];
        }

        private void calcktrace(int iter_num, int discrete_num)
        {
            double sum = 0;
            for (int i = 0; i < this.system_size; i++)
                sum += this.A_matrix[iter_num, discrete_num, i, i];
            this.P_trace[iter_num, discrete_num] = sum / (iter_num + 1);
        }

    }
}
